6  PBMCs: Determine IMP Expanded Clones

6.1 Set up workspace

# Libraries
library(dplyr)

Attaching package: 'dplyr'

The following objects are masked from 'package:stats':

    filter, lag

The following objects are masked from 'package:base':

    intersect, setdiff, setequal, union

library(tidyverse)

── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
✔ forcats   1.0.0     ✔ readr     2.1.5
✔ ggplot2   3.5.1     ✔ stringr   1.5.1
✔ lubridate 1.9.4     ✔ tibble    3.2.1
✔ purrr     1.0.4     ✔ tidyr     1.3.1

── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag()    masks stats::lag()
ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors

library(ggplot2)
library(patchwork)
library(fossil)

Loading required package: sp
Loading required package: maps

Attaching package: 'maps'

The following object is masked from 'package:purrr':

    map

Loading required package: shapefiles
Loading required package: foreign

Attaching package: 'shapefiles'

The following objects are masked from 'package:foreign':

    read.dbf, write.dbf

library(ggridges)
library(webr)
library(scales)

Attaching package: 'scales'

The following object is masked from 'package:purrr':

    discard

The following object is masked from 'package:readr':

    col_factor

# Function to calculate lower limit of detection
lod <- function(dat){
  10^(floor(log10(min(dat, na.rm = TRUE))))
}

# Function to calculate the minimum value in a bimodal distribution
find_minimum <- function(fc_vector, adjust = 1){
  # Write density function
  dx <- density(fc_vector, adjust = adjust)
  # Calculate first derivative
  slope <- (diff(dx$y)/diff(dx$x))
  # Find where the derivatives change sign
  change_rows <- which(diff(sign(slope))!=0)
  change_x <- dx$x[change_rows]
  # Extract the second time it changes sign
  min_x <- change_x[2]
  # second_min_x <- change_x[4]
  return(min_x)
}

6.2 Categorizations of temporally dynamic clones

Defintion of Post-Vaccine clones in PBMCs 1. Quality control metric: the clone must have 3 or more UMI in at least one time point 2. Below the limit of detection pre-treatment and Pre-vaccine 3. >10x lower limit at post-vaccine 4. Present at all time points post-vaccine

Defintion of Post-Nivolumab expanded clones in PBMCs 1. Quality control metric: the clone must have 3 or more UMI in at least one time point 2. Below the limit of detection post-nivolumab 3. >10x lower limit post-nivolumab 4. Present at all time points after post-nivolumab

Defintion of Existing clones in PBMCs 1. Quality control metric: the clone must have 3 or more UMI in at least one time point 2. >10x lower limit pre-treatment 3. Present at all time points

6.3 Load PBMC data

p101_betas <- read.csv("p101_betas_merged_Part1.csv")
p103_betas <- read.csv("p103_betas_merged_Part1.csv")
p104_betas <- read.csv("p104_betas_merged_Part1.csv")
p105_betas <- read.csv("p105_betas_merged_Part1.csv")
p106_betas <- read.csv("p106_betas_merged_Part1.csv")
p108_betas <- read.csv("p108_betas_merged_Part1.csv")
p109_betas <- read.csv("p109_betas_merged_Part1.csv")
p110_betas <- read.csv("p110_betas_merged_Part1.csv")
p111_betas <- read.csv("p111_betas_merged_Part1.csv")

6.4 Change all undetected longitudinal values to the lower limit of detection per timepoint

Makes plotting the movement of clones over time easier!

# Generate the lower limit of detection before any filters are applied
p101_pretreatment_lod <- lod(p101_betas$p101_pretreatment)
p101_prevax_lod <- lod(p101_betas$p101_prevax)
p101_postvax_lod <- lod(p101_betas$p101_postvax)

p103_pretreatment_lod <- lod(p103_betas$p103_pretreatment)
p103_prevax_lod <- lod(p103_betas$p103_prevax)
p103_postvax_lod <- lod(p103_betas$p103_postvax)
p103_w48_lod <- lod(p103_betas$p103_w48)
p103_w72_lod <- lod(p103_betas$p103_w72)

p104_pretreatment_lod <- lod(p104_betas$p104_pretreatment)
p104_prevax_lod <- lod(p104_betas$p104_prevax)
p104_postvax_lod <- lod(p104_betas$p104_postvax)
p104_w48_lod <- lod(p104_betas$p104_w48)

p105_pretreatment_lod <- lod(p105_betas$p105_pretreatment)
p105_prevax_lod <- lod(p105_betas$p105_prevax)
p105_postvax_lod <- lod(p105_betas$p105_postvax)
p105_w48_lod <- lod(p105_betas$p105_w48)

p106_pretreatment_lod <- lod(p106_betas$p106_pretreatment)
p106_prevax_lod <- lod(p106_betas$p106_prevax)
p106_postvax_lod <- lod(p106_betas$p106_postvax)
p106_w48_lod <- lod(p106_betas$p106_w48)

p108_pretreatment_lod <- lod(p108_betas$p108_pretreatment)
p108_prevax_lod <- lod(p108_betas$p108_prevax)
p108_postvax_lod <- lod(p108_betas$p108_postvax)
p108_w32_lod <- lod(p108_betas$p108_w32)

p109_pretreatment_lod <- lod(p109_betas$p109_pretreatment)
p109_prevax_lod <- lod(p109_betas$p109_prevax)
p109_postvax_lod <- lod(p109_betas$p109_postvax)
p109_w48_lod <- lod(p109_betas$p109_w48)

p110_pretreatment_lod <- lod(p110_betas$p110_pretreatment)
p110_prevax_lod <- lod(p110_betas$p110_prevax)
p110_postvax_lod <- lod(p110_betas$p110_postvax)
p110_w38_lod <- lod(p110_betas$p110_w38)

p111_pretreatment_lod <- lod(p111_betas$p111_pretreatment)
p111_prevax_lod <- lod(p111_betas$p111_prevax)
p111_postvax_lod <- lod(p111_betas$p111_postvax)
p111_w40_lod <- lod(p111_betas$p111_w40)

# Update undetected values to the LOD
p101_betas <- p101_betas %>%
  mutate(p101_pretreatment = replace_na(p101_pretreatment, p101_pretreatment_lod),
         p101_prevax = replace_na(p101_prevax, p101_prevax_lod),
         p101_postvax = replace_na(p101_postvax, p101_postvax_lod))

p103_betas <- p103_betas %>%
  mutate(p103_pretreatment = replace_na(p103_pretreatment, p103_pretreatment_lod),
         p103_prevax = replace_na(p103_prevax, p103_prevax_lod),
         p103_postvax = replace_na(p103_postvax, p103_postvax_lod),
         p103_w48 = replace_na(p103_w48, p103_w48_lod),
         p103_w72 = replace_na(p103_w72, p103_w72_lod))
p104_betas <- p104_betas %>%
  mutate(p104_pretreatment = replace_na(p104_pretreatment, p104_pretreatment_lod),
         p104_prevax = replace_na(p104_prevax, p104_prevax_lod),
         p104_postvax = replace_na(p104_postvax, p104_postvax_lod),
         p104_w48 = replace_na(p104_w48, p104_w48_lod))
p105_betas <- p105_betas %>%
  mutate(p105_pretreatment = replace_na(p105_pretreatment, p105_pretreatment_lod),
         p105_prevax = replace_na(p105_prevax, p105_prevax_lod),
         p105_postvax = replace_na(p105_postvax, p105_postvax_lod),
         p105_w48 = replace_na(p105_w48, p105_w48_lod))
p106_betas <- p106_betas %>%
  mutate(p106_pretreatment = replace_na(p106_pretreatment, p106_pretreatment_lod),
         p106_prevax = replace_na(p106_prevax, p106_prevax_lod),
         p106_postvax = replace_na(p106_postvax, p106_postvax_lod),
         p106_w48 = replace_na(p106_w48, p106_w48_lod))
p108_betas <- p108_betas %>%
  mutate(p108_pretreatment = replace_na(p108_pretreatment, p108_pretreatment_lod),
         p108_prevax = replace_na(p108_prevax, p108_prevax_lod),
         p108_postvax = replace_na(p108_postvax, p108_postvax_lod),
         p108_w32 = replace_na(p108_w32, p108_w32_lod))
p109_betas <- p109_betas %>%
  mutate(p109_pretreatment = replace_na(p109_pretreatment, p109_pretreatment_lod),
         p109_prevax = replace_na(p109_prevax, p109_prevax_lod),
         p109_postvax = replace_na(p109_postvax, p109_postvax_lod),
         p109_w48 = replace_na(p109_w48, p109_w48_lod))
p110_betas <- p110_betas %>%
  mutate(p110_pretreatment = replace_na(p110_pretreatment, p110_pretreatment_lod),
         p110_prevax = replace_na(p110_prevax, p110_prevax_lod),
         p110_postvax = replace_na(p110_postvax, p110_postvax_lod),
         p110_w38 = replace_na(p110_w38, p110_w38_lod))
p111_betas <- p111_betas %>%
  mutate(p111_pretreatment = replace_na(p111_pretreatment, p111_pretreatment_lod),
         p111_prevax = replace_na(p111_prevax, p111_prevax_lod),
         p111_postvax = replace_na(p111_postvax, p111_postvax_lod),
         p111_w40 = replace_na(p111_w40, p111_w40_lod))

6.5 Load categorized clones

preexisting_betas <- read.csv("PBMC_Existing_Beta_Chains_Part2.csv", row.names = 1)
postnivo_betas <- read.csv("PBMC_Post-Nivolumab_Beta_Chains_Part2.csv", row.names = 1)
postvax_betas <- read.csv("PBMC_Post-Vaccine_Beta_Chains_Part2.csv", row.names = 1)

# Check how many beta chains are shared across patients
table(duplicated(postvax_betas$Beta_clonotype))

FALSE  TRUE 
 8656    25 

table(duplicated(postnivo_betas$Beta_clonotype))

FALSE  TRUE 
 5621    12 

table(duplicated(preexisting_betas$Beta_clonotype))

FALSE  TRUE 
49282   337 

6.6 Integrate categories into data, adding an additional category called “Other” for clones who do not match the other criteria

p101_betas_typed <- p101_betas %>%
  mutate(Type = case_when(Beta_clonotype %in% filter(postvax_betas, postvax_betas$Patient == "P101")$Beta_clonotype ~ "Post-Vaccine",
                          Beta_clonotype %in% filter(postnivo_betas, postnivo_betas$Patient == "P101")$Beta_clonotype ~ "Post-Nivolumab",
                          Beta_clonotype %in% filter(preexisting_betas, preexisting_betas$Patient == "P101")$Beta_clonotype ~ "Existing"),
         Type = case_when(!is.na(Type) ~ Type,
                          is.na(Type) ~ "Other"),
         Patient = "P101")

p103_betas_typed <- p103_betas %>%
  mutate(Type = case_when(Beta_clonotype %in% filter(postvax_betas, postvax_betas$Patient == "P103")$Beta_clonotype ~ "Post-Vaccine",
                          Beta_clonotype %in% filter(postnivo_betas, postnivo_betas$Patient == "P103")$Beta_clonotype ~ "Post-Nivolumab",
                          Beta_clonotype %in% filter(preexisting_betas, preexisting_betas$Patient == "P103")$Beta_clonotype ~ "Existing"),
         Type = case_when(!is.na(Type) ~ Type,
                          is.na(Type) ~ "Other"),
         Patient = "P103")

p104_betas_typed <- p104_betas %>%
  mutate(Type = case_when(Beta_clonotype %in% filter(postvax_betas, postvax_betas$Patient == "P104")$Beta_clonotype ~ "Post-Vaccine",
                          Beta_clonotype %in% filter(postnivo_betas, postnivo_betas$Patient == "P104")$Beta_clonotype ~ "Post-Nivolumab",
                          Beta_clonotype %in% filter(preexisting_betas, preexisting_betas$Patient == "P104")$Beta_clonotype ~ "Existing"),
         Type = case_when(!is.na(Type) ~ Type,
                          is.na(Type) ~ "Other"),
         Patient = "P104")

p105_betas_typed <- p105_betas %>%
  mutate(Type = case_when(Beta_clonotype %in% filter(postvax_betas, postvax_betas$Patient == "P105")$Beta_clonotype ~ "Post-Vaccine",
                          Beta_clonotype %in% filter(postnivo_betas, postnivo_betas$Patient == "P105")$Beta_clonotype ~ "Post-Nivolumab",
                          Beta_clonotype %in% filter(preexisting_betas, preexisting_betas$Patient == "P105")$Beta_clonotype ~ "Existing"),
         Type = case_when(!is.na(Type) ~ Type,
                          is.na(Type) ~ "Other"),
         Patient = "P105")

p106_betas_typed <- p106_betas %>%
  mutate(Type = case_when(Beta_clonotype %in% filter(postvax_betas, postvax_betas$Patient == "P106")$Beta_clonotype ~ "Post-Vaccine",
                          Beta_clonotype %in% filter(postnivo_betas, postnivo_betas$Patient == "P106")$Beta_clonotype ~ "Post-Nivolumab",
                          Beta_clonotype %in% filter(preexisting_betas, preexisting_betas$Patient == "P106")$Beta_clonotype ~ "Existing"),
         Type = case_when(!is.na(Type) ~ Type,
                          is.na(Type) ~ "Other"),
         Patient = "P106")

p108_betas_typed <- p108_betas %>%
  mutate(Type = case_when(Beta_clonotype %in% filter(postvax_betas, postvax_betas$Patient == "P108")$Beta_clonotype ~ "Post-Vaccine",
                          Beta_clonotype %in% filter(postnivo_betas, postnivo_betas$Patient == "P108")$Beta_clonotype ~ "Post-Nivolumab",
                          Beta_clonotype %in% filter(preexisting_betas, preexisting_betas$Patient == "P108")$Beta_clonotype ~ "Existing"),
         Type = case_when(!is.na(Type) ~ Type,
                          is.na(Type) ~ "Other"),
         Patient = "P108")

p109_betas_typed <- p109_betas %>%
  mutate(Type = case_when(Beta_clonotype %in% filter(postvax_betas, postvax_betas$Patient == "P109")$Beta_clonotype ~ "Post-Vaccine",
                          Beta_clonotype %in% filter(postnivo_betas, postnivo_betas$Patient == "P109")$Beta_clonotype ~ "Post-Nivolumab",
                          Beta_clonotype %in% filter(preexisting_betas, preexisting_betas$Patient == "P109")$Beta_clonotype ~ "Existing"),
         Type = case_when(!is.na(Type) ~ Type,
                          is.na(Type) ~ "Other"),
         Patient = "P109")

p110_betas_typed <- p110_betas %>%
  mutate(Type = case_when(Beta_clonotype %in% filter(postvax_betas, postvax_betas$Patient == "P110")$Beta_clonotype ~ "Post-Vaccine",
                          Beta_clonotype %in% filter(postnivo_betas, postnivo_betas$Patient == "P110")$Beta_clonotype ~ "Post-Nivolumab",
                          Beta_clonotype %in% filter(preexisting_betas, preexisting_betas$Patient == "P110")$Beta_clonotype ~ "Existing"),
         Type = case_when(!is.na(Type) ~ Type,
                          is.na(Type) ~ "Other"),
         Patient = "P110")

p111_betas_typed <- p111_betas %>%
  mutate(Type = case_when(Beta_clonotype %in% filter(postvax_betas, postvax_betas$Patient == "P111")$Beta_clonotype ~ "Post-Vaccine",
                          Beta_clonotype %in% filter(postnivo_betas, postnivo_betas$Patient == "P111")$Beta_clonotype ~ "Post-Nivolumab",
                          Beta_clonotype %in% filter(preexisting_betas, preexisting_betas$Patient == "P111")$Beta_clonotype ~ "Existing"),
         Type = case_when(!is.na(Type) ~ Type,
                          is.na(Type) ~ "Other"),
         Patient = "P111")

6.7 Save categorized clones

write.csv(p101_betas_typed, "p101_betas_merged_typed_full_Part3.csv", row.names = FALSE)
write.csv(p103_betas_typed, "p103_betas_merged_typed_full_Part3.csv", row.names = FALSE)
write.csv(p104_betas_typed, "p104_betas_merged_typed_full_Part3.csv", row.names = FALSE)
write.csv(p105_betas_typed, "p105_betas_merged_typed_full_Part3.csv", row.names = FALSE)
write.csv(p106_betas_typed, "p106_betas_merged_typed_full_Part3.csv", row.names = FALSE)
write.csv(p108_betas_typed, "p108_betas_merged_typed_full_Part3.csv", row.names = FALSE)
write.csv(p109_betas_typed, "p109_betas_merged_typed_full_Part3.csv", row.names = FALSE)
write.csv(p110_betas_typed, "p110_betas_merged_typed_full_Part3.csv", row.names = FALSE)
write.csv(p111_betas_typed, "p111_betas_merged_typed_full_Part3.csv", row.names = FALSE)

6.8 Remove clones deemed “Other” from data

p101_betas_typed <- p101_betas_typed %>%
  filter(Type != "Other") %>%
  mutate(Type = factor(Type, levels = c("Post-Vaccine", "Post-Nivolumab", "Existing")))
p103_betas_typed <- p103_betas_typed %>%
  filter(Type != "Other") %>%
  mutate(Type = factor(Type, levels = c("Post-Vaccine", "Post-Nivolumab", "Existing")))
p104_betas_typed <- p104_betas_typed %>%
  filter(Type != "Other") %>%
  mutate(Type = factor(Type, levels = c("Post-Vaccine", "Post-Nivolumab", "Existing")))
p105_betas_typed <- p105_betas_typed %>%
  filter(Type != "Other") %>%
  mutate(Type = factor(Type, levels = c("Post-Vaccine", "Post-Nivolumab", "Existing")))
p106_betas_typed <- p106_betas_typed %>%
  filter(Type != "Other") %>%
  mutate(Type = factor(Type, levels = c("Post-Vaccine", "Post-Nivolumab", "Existing")))
p108_betas_typed <- p108_betas_typed %>%
  filter(Type != "Other") %>%
  mutate(Type = factor(Type, levels = c("Post-Vaccine", "Post-Nivolumab", "Existing")))
p109_betas_typed <- p109_betas_typed %>%
  filter(Type != "Other") %>%
  mutate(Type = factor(Type, levels = c("Post-Vaccine", "Post-Nivolumab", "Existing")))
p110_betas_typed <- p110_betas_typed %>%
  filter(Type != "Other") %>%
  mutate(Type = factor(Type, levels = c("Post-Vaccine", "Post-Nivolumab", "Existing")))
p111_betas_typed <- p111_betas_typed %>%
  filter(Type != "Other") %>%
  mutate(Type = factor(Type, levels = c("Post-Vaccine", "Post-Nivolumab", "Existing")))

6.9 Plot ridgeplot split by patient and TCR type

p101_ridge <- p101_betas_typed %>%
  ggplot(aes(x = log2(max_fc_pool_vs_dmso), y = Type, fill = Type)) +
  geom_density_ridges(bandwidth = 0.3, scale = 2) +
  # coord_cartesian(xlim = c(1.5, 1000)) +
  scale_fill_manual(values = c("#0072B2", "darkgoldenrod3", "grey70")) +
  # scale_x_log10(breaks=c(1.5, 10, 100, 1000)) +
  xlab("Invitro Expansion") +
  ylab("Number of clones") +
  ggtitle("P101") +
  theme_classic() +
  theme(plot.title = element_text(hjust=0.5))

p103_ridge <- p103_betas_typed %>%
  ggplot(aes(x = log2(max_fc_pool_vs_dmso), y = Type, fill = Type)) +
  geom_density_ridges(bandwidth = 0.3, scale = 2) +
  # coord_cartesian(xlim = c(1.5, 1000)) +
  scale_fill_manual(values = c("#0072B2", "darkgoldenrod3", "grey70")) +
  # scale_x_log10(breaks=c(1.5, 10, 100, 1000)) +
  xlab("Invitro Expansion") +
  ylab("Number of clones") +
  ggtitle("P103") +
  theme_classic() +
  theme(plot.title = element_text(hjust=0.5))

p104_ridge <- p104_betas_typed %>%
  ggplot(aes(x = log2(max_fc_pool_vs_dmso), y = Type, fill = Type)) +
  geom_density_ridges(bandwidth = 0.3, scale = 2) +
  # coord_cartesian(xlim = c(1.5, 1000)) +
  scale_fill_manual(values = c("#0072B2", "darkgoldenrod3", "grey70")) +
  # scale_x_log10(breaks=c(1.5, 10, 100, 1000)) +
  xlab("Invitro Expansion") +
  ylab("Number of clones") +
  ggtitle("P104") +
  theme_classic() +
  theme(plot.title = element_text(hjust=0.5))

p105_ridge <- p105_betas_typed %>%
  ggplot(aes(x = log2(max_fc_pool_vs_dmso), y = Type, fill = Type)) +
  geom_density_ridges(bandwidth = 0.3, scale = 2) +
  # coord_cartesian(xlim = c(1.5, 1000)) +
  scale_fill_manual(values = c("#0072B2", "darkgoldenrod3", "grey70")) +
  # scale_x_log10(breaks=c(1.5, 10, 100, 1000)) +
  xlab("Invitro Expansion") +
  ylab("Number of clones") +
  ggtitle("P105") +
  theme_classic() +
  theme(plot.title = element_text(hjust=0.5))

p106_ridge <- p106_betas_typed %>%
  ggplot(aes(x = log2(max_fc_pool_vs_dmso), y = Type, fill = Type)) +
  geom_density_ridges(bandwidth = 0.3, scale = 2) +
  # coord_cartesian(xlim = c(1.5, 1000)) +
  scale_fill_manual(values = c("#0072B2", "darkgoldenrod3", "grey70")) +
  # scale_x_log10(breaks=c(1.5, 10, 100, 1000)) +
  xlab("Invitro Expansion") +
  ylab("Number of clones") +
  ggtitle("P106") +
  theme_classic() +
  theme(plot.title = element_text(hjust=0.5))

p108_ridge <- p108_betas_typed %>%
  ggplot(aes(x = log2(max_fc_pool_vs_dmso), y = Type, fill = Type)) +
  geom_density_ridges(bandwidth = 0.3, scale = 2) +
  # coord_cartesian(xlim = c(1.5, 1000)) +
  scale_fill_manual(values = c("#0072B2", "darkgoldenrod3", "grey70")) +
  # scale_x_log10(breaks=c(1.5, 10, 100, 1000)) +
  xlab("Invitro Expansion") +
  ylab("Number of clones") +
  ggtitle("P108") +
  theme_classic() +
  theme(plot.title = element_text(hjust=0.5))

p109_ridge <- p109_betas_typed %>%
  ggplot(aes(x = log2(max_fc_pool_vs_dmso), y = Type, fill = Type)) +
  geom_density_ridges(bandwidth = 0.3, scale = 2) +
  # coord_cartesian(xlim = c(1.5, 1000)) +
  scale_fill_manual(values = c("#0072B2", "darkgoldenrod3", "grey70")) +
  # scale_x_log10(breaks=c(1.5, 10, 100, 1000)) +
  xlab("Invitro Expansion") +
  ylab("Number of clones") +
  ggtitle("P109") +
  theme_classic() +
  theme(plot.title = element_text(hjust=0.5))

p110_ridge <- p110_betas_typed %>%
  ggplot(aes(x = log2(max_fc_pool_vs_dmso), y = Type, fill = Type)) +
  geom_density_ridges(bandwidth = 0.3, scale = 2) +
  # coord_cartesian(xlim = c(1.5, 1000)) +
  scale_fill_manual(values = c("#0072B2", "darkgoldenrod3", "grey70")) +
  # scale_x_log10(breaks=c(1.5, 10, 100, 1000)) +
  xlab("Invitro Expansion") +
  ylab("Number of clones") +
  ggtitle("P110") +
  theme_classic() +
  theme(plot.title = element_text(hjust=0.5))

p111_ridge <- p111_betas_typed %>%
  ggplot(aes(x = log2(max_fc_pool_vs_dmso), y = Type, fill = Type)) +
  geom_density_ridges(bandwidth = 0.3, scale = 2) +
  # coord_cartesian(xlim = c(1.5, 1000)) +
  scale_fill_manual(values = c("#0072B2", "darkgoldenrod3", "grey70")) +
  # scale_x_log10(breaks=c(1.5, 10, 100, 1000)) +
  xlab("Invitro Expansion") +
  ylab("Number of clones") +
  ggtitle("P111") +
  theme_classic() +
  theme(plot.title = element_text(hjust=0.5))

p101_ridge + p103_ridge + p104_ridge + p105_ridge + p106_ridge +
  p108_ridge + p109_ridge + p110_ridge + p111_ridge + plot_layout(guides = 'collect', 
                                                                  axis_titles = "collect", 
                                                                  ncol = 3)

6.10 Find a cutoff for invitro expansion per patient

p101_max_fc <- p101_betas_typed$max_fc_pool_vs_dmso[!is.na(p101_betas_typed$max_fc_pool_vs_dmso)]
p103_max_fc <- p103_betas_typed$max_fc_pool_vs_dmso[!is.na(p103_betas_typed$max_fc_pool_vs_dmso)]
p104_max_fc <- p104_betas_typed$max_fc_pool_vs_dmso[!is.na(p104_betas_typed$max_fc_pool_vs_dmso)]
p105_max_fc <- p105_betas_typed$max_fc_pool_vs_dmso[!is.na(p105_betas_typed$max_fc_pool_vs_dmso)]
p106_max_fc <- p106_betas_typed$max_fc_pool_vs_dmso[!is.na(p106_betas_typed$max_fc_pool_vs_dmso)]
p108_max_fc <- p108_betas_typed$max_fc_pool_vs_dmso[!is.na(p108_betas_typed$max_fc_pool_vs_dmso)]
p109_max_fc <- p109_betas_typed$max_fc_pool_vs_dmso[!is.na(p109_betas_typed$max_fc_pool_vs_dmso)]
p110_max_fc <- p110_betas_typed$max_fc_pool_vs_dmso[!is.na(p110_betas_typed$max_fc_pool_vs_dmso)]
p111_max_fc <- p111_betas_typed$max_fc_pool_vs_dmso[!is.na(p111_betas_typed$max_fc_pool_vs_dmso)]

p101_min_x <- find_minimum(p101_max_fc)
p103_min_x <- find_minimum(p103_max_fc)
p104_min_x <- find_minimum(p104_max_fc)
p105_min_x <- find_minimum(p105_max_fc)
p106_min_x <- find_minimum(p106_max_fc)
p108_min_x <- find_minimum(p108_max_fc)
p109_min_x <- find_minimum(p109_max_fc)
p110_min_x <- find_minimum(p110_max_fc)
p111_min_x <- find_minimum(p111_max_fc)

6.11 Ridgeplot with expansion minimum

all_pt_ridge <- (p101_ridge + geom_vline(xintercept = log2(p101_min_x), color = "red", linetype="dashed")) + 
  (p103_ridge + geom_vline(xintercept = log2(p103_min_x), color = "red", linetype="dashed")) + 
  (p104_ridge + geom_vline(xintercept = log2(p104_min_x), color = "red", linetype="dashed")) + 
  (p105_ridge + geom_vline(xintercept = log2(p105_min_x), color = "red", linetype="dashed")) + 
  (p106_ridge + geom_vline(xintercept = log2(p106_min_x), color = "red", linetype="dashed")) +
  (p108_ridge + geom_vline(xintercept = log2(p108_min_x), color = "red", linetype="dashed")) + 
  (p109_ridge + geom_vline(xintercept = log2(p109_min_x), color = "red", linetype="dashed")) + 
  (p110_ridge + geom_vline(xintercept = log2(p110_min_x), color = "red", linetype="dashed")) + 
  (p111_ridge + geom_vline(xintercept = log2(p111_min_x), color = "red", linetype="dashed")) + 
  plot_layout(guides = 'collect', axis_titles = "collect", ncol = 3)

all_pt_ridge

6.12 Which clones see invitro expansion past local minimum?

p101_betas_typed <- p101_betas_typed %>%
  mutate(expanded_local_min = case_when(max_fc_pool_vs_dmso >= p101_min_x ~ TRUE,
                                        max_fc_pool_vs_dmso < p101_min_x ~ FALSE),
         Patient = "P101")
p103_betas_typed <- p103_betas_typed %>%
  mutate(expanded_local_min = case_when(max_fc_pool_vs_dmso >= p103_min_x ~ TRUE,
                                        max_fc_pool_vs_dmso < p103_min_x ~ FALSE),
         Patient = "P103")
p104_betas_typed <- p104_betas_typed %>%
  mutate(expanded_local_min = case_when(max_fc_pool_vs_dmso >= p104_min_x ~ TRUE,
                                        max_fc_pool_vs_dmso < p104_min_x ~ FALSE),
         Patient = "P104")
p105_betas_typed <- p105_betas_typed %>%
  mutate(expanded_local_min = case_when(max_fc_pool_vs_dmso >= p105_min_x ~ TRUE,
                                        max_fc_pool_vs_dmso < p105_min_x ~ FALSE),
         Patient = "P105")
p106_betas_typed <- p106_betas_typed %>%
  mutate(expanded_local_min = case_when(max_fc_pool_vs_dmso >= p106_min_x ~ TRUE,
                                        max_fc_pool_vs_dmso < p106_min_x ~ FALSE),
         Patient = "P106")
p108_betas_typed <- p108_betas_typed %>%
  mutate(expanded_local_min = case_when(max_fc_pool_vs_dmso >= p108_min_x ~ TRUE,
                                        max_fc_pool_vs_dmso < p108_min_x ~ FALSE),
         Patient = "P108")
p109_betas_typed <- p109_betas_typed %>%
  mutate(expanded_local_min = case_when(max_fc_pool_vs_dmso >= p109_min_x ~ TRUE,
                                        max_fc_pool_vs_dmso < p109_min_x ~ FALSE),
         Patient = "P109")
p110_betas_typed <- p110_betas_typed %>%
  mutate(expanded_local_min = case_when(max_fc_pool_vs_dmso >= p110_min_x ~ TRUE,
                                        max_fc_pool_vs_dmso < p110_min_x ~ FALSE),
         Patient = "P110")
p111_betas_typed <- p111_betas_typed %>%
  mutate(expanded_local_min = case_when(max_fc_pool_vs_dmso >= p111_min_x ~ TRUE,
                                        max_fc_pool_vs_dmso < p111_min_x ~ FALSE),
         Patient = "P111")

6.13 Subset for clones which pass invitro expansion- i.e. are “IMP expanded clones”

p101_betas_imp_expanded <- p101_betas_typed %>%
  filter(expanded_local_min == TRUE)
p103_betas_imp_expanded <- p103_betas_typed %>%
  filter(expanded_local_min == TRUE)
p104_betas_imp_expanded <- p104_betas_typed %>%
  filter(expanded_local_min == TRUE)
p105_betas_imp_expanded <- p105_betas_typed %>%
  filter(expanded_local_min == TRUE)
p106_betas_imp_expanded <- p106_betas_typed %>%
  filter(expanded_local_min == TRUE)
p108_betas_imp_expanded <- p108_betas_typed %>%
  filter(expanded_local_min == TRUE)
p109_betas_imp_expanded <- p109_betas_typed %>%
  filter(expanded_local_min == TRUE)
p110_betas_imp_expanded <- p110_betas_typed %>%
  filter(expanded_local_min == TRUE)
p111_betas_imp_expanded <- p111_betas_typed %>%
  filter(expanded_local_min == TRUE)

6.14 Save typed clones

write.csv(p101_betas_typed, "p101_betas_merged_typed_Part3.csv", row.names = FALSE)
write.csv(p103_betas_typed, "p103_betas_merged_typed_Part3.csv", row.names = FALSE)
write.csv(p104_betas_typed, "p104_betas_merged_typed_Part3.csv", row.names = FALSE)
write.csv(p105_betas_typed, "p105_betas_merged_typed_Part3.csv", row.names = FALSE)
write.csv(p106_betas_typed, "p106_betas_merged_typed_Part3.csv", row.names = FALSE)
write.csv(p108_betas_typed, "p108_betas_merged_typed_Part3.csv", row.names = FALSE)
write.csv(p109_betas_typed, "p109_betas_merged_typed_Part3.csv", row.names = FALSE)
write.csv(p110_betas_typed, "p110_betas_merged_typed_Part3.csv", row.names = FALSE)
write.csv(p111_betas_typed, "p111_betas_merged_typed_Part3.csv", row.names = FALSE)

6.15 Save expanded clones

write.csv(p101_betas_imp_expanded, "p101_betas_merged_typed_imp_expanded_Part3.csv", row.names = FALSE)
write.csv(p103_betas_imp_expanded, "p103_betas_merged_typed_imp_expanded_Part3.csv", row.names = FALSE)
write.csv(p104_betas_imp_expanded, "p104_betas_merged_typed_imp_expanded_Part3.csv", row.names = FALSE)
write.csv(p105_betas_imp_expanded, "p105_betas_merged_typed_imp_expanded_Part3.csv", row.names = FALSE)
write.csv(p106_betas_imp_expanded, "p106_betas_merged_typed_imp_expanded_Part3.csv", row.names = FALSE)
write.csv(p108_betas_imp_expanded, "p108_betas_merged_typed_imp_expanded_Part3.csv", row.names = FALSE)
write.csv(p109_betas_imp_expanded, "p109_betas_merged_typed_imp_expanded_Part3.csv", row.names = FALSE)
write.csv(p110_betas_imp_expanded, "p110_betas_merged_typed_imp_expanded_Part3.csv", row.names = FALSE)
write.csv(p111_betas_imp_expanded, "p111_betas_merged_typed_imp_expanded_Part3.csv", row.names = FALSE)

6.16 Save the invitro-expansion cutoff values per patient

pt_vec <- c("P101", "P103", "P104", "P105", "P106", "P108", "P109", "P110", "P111")
cutoff_vec <- c(p101_min_x, p103_min_x, p104_min_x, p105_min_x, p106_min_x, p108_min_x, p109_min_x, p110_min_x, p111_min_x)
cutoff_df <- data.frame(Patient = pt_vec, invitro_exp_cutoff = cutoff_vec)

write.csv(cutoff_df, "IMP_expanded_cuotff_Part3.csv", row.names = FALSE)

6.17 Get session info

sessionInfo()

R version 4.3.2 (2023-10-31)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Rocky Linux 8.10 (Green Obsidian)

Matrix products: default
BLAS/LAPACK: /usr/lib64/libopenblasp-r0.3.15.so;  LAPACK version 3.9.0

locale:
 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
 [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
 [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
 [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C            
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       

time zone: America/New_York
tzcode source: system (glibc)

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] scales_1.3.0     webr_0.1.5       ggridges_0.5.6   fossil_0.4.0    
 [5] shapefiles_0.7.2 foreign_0.8-90   maps_3.4.2.1     sp_2.2-0        
 [9] patchwork_1.3.0  lubridate_1.9.4  forcats_1.0.0    stringr_1.5.1   
[13] purrr_1.0.4      readr_2.1.5      tidyr_1.3.1      tibble_3.2.1    
[17] ggplot2_3.5.1    tidyverse_2.0.0  dplyr_1.1.4     

loaded via a namespace (and not attached):
 [1] mnormt_2.1.1            rlang_1.1.5             magrittr_2.0.3         
 [4] rio_1.2.3               compiler_4.3.2          systemfonts_1.2.2      
 [7] vctrs_0.6.5             devEMF_4.5              pkgconfig_2.0.3        
[10] crayon_1.5.3            fastmap_1.2.0           backports_1.5.0        
[13] labeling_0.4.3          promises_1.3.2          rmarkdown_2.29         
[16] tzdb_0.5.0              ragg_1.3.3              xfun_0.50              
[19] cachem_1.1.0            jsonlite_1.8.9          later_1.4.1            
[22] sjmisc_2.8.10           uuid_1.2-1              tweenr_2.0.3           
[25] psych_2.5.3             broom_1.0.8             parallel_4.3.2         
[28] R6_2.6.1                vcd_1.4-13              RColorBrewer_1.1-3     
[31] bslib_0.8.0             stringi_1.8.4           car_3.1-3              
[34] lmtest_0.9-40           jquerylib_0.1.4         Rcpp_1.0.14            
[37] knitr_1.49              zoo_1.8-13              httpuv_1.6.15          
[40] timechange_0.3.0        tidyselect_1.2.1        rstudioapi_0.17.1      
[43] abind_1.4-8             moonBook_0.3.1          miniUI_0.1.1.1         
[46] sjlabelled_1.2.0        lattice_0.22-7          shiny_1.9.1            
[49] rrtable_0.3.0           withr_3.0.2             flextable_0.9.7        
[52] askpass_1.2.1           evaluate_1.0.1          polyclip_1.10-7        
[55] zip_2.3.2               xml2_1.3.8              pillar_1.10.1          
[58] ggpubr_0.6.0            carData_3.0-5           DT_0.33                
[61] insight_0.20.3          generics_0.1.3          hms_1.1.3              
[64] munsell_0.5.1           xtable_1.8-4            glue_1.8.0             
[67] gdtools_0.4.2           tools_4.3.2             data.table_1.15.4      
[70] openxlsx_4.2.8          ggsignif_0.6.4          grid_4.3.2             
[73] editData_0.1.8          colorspace_2.1-1        nlme_3.1-168           
[76] ggforce_0.4.2           Formula_1.2-5           cli_3.6.3              
[79] textshaping_1.0.0       rvg_0.3.3               officer_0.6.8          
[82] fontBitstreamVera_0.1.1 ztable_0.2.3            gtable_0.3.6           
[85] rstatix_0.7.2           sass_0.4.9              digest_0.6.37          
[88] fontquiver_0.2.1        htmlwidgets_1.6.4       farver_2.1.2           
[91] htmltools_0.5.8.1       lifecycle_1.0.4         shinyWidgets_0.9.0     
[94] mime_0.13               fontLiberation_0.1.0    openssl_2.3.1          
[97] MASS_7.3-60.0.1